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ABSTRACT 



The development of a numerical 1000-mb prognosis by means 
of a 500 to 1000-mb thickness forecast is investigated. The 
thickness forecast is based, in part, on the 500 -mb barotropic 
prognosis as developed by Fleet Numerical Weather Facility, 
Monterey. 

Forecasts were prepared for three days during May 1966, 
with various amounts of smoothing of the component fields. 

The model showed promise of producing a rapid, usable prognosis; 
however further testing is needed with the probable addition of 
some climatic constraints. It is also hoped that the yet untried 
diabatic heating and surface friction terms will provide some 



further improvement 
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1 . Introduction 



With the development of the high speed digital computer, the 
problem of weather prognosis has evolved from one of qualitative 
reasoning by an experienced forecaster to one of objective 
computations based on mathematical equations representing the 
motions of the atmosphere. At this time, daily numerical analysis 
and prognostic charts are being produced at Fleet Numerical Weather 
Facility, Monterey (FNWF) , and rapidly disseminated to fleet units 
around the world. The entire field of meteorology is in a rapid state 
of change with more and more products being produced numerically. 

The most successful numerical prognosis has been at the 500 mb 
level. If a successful numerical isobaric thickness forecast could 
be developed, the success of the 500-mb forecast could then be 
reflected to higher and lower levels. The objective of this research 
is to obtain a dependable 24-hour numerical forecast of thickness 
for the layer between 500-mb and 1000-mb, and then using FNWF's 
500-mb barotropic forecast obtain a 1000-mb height prognosis. 

The thickness of an air layer at any given time and place is 
dependent upon the mean temperature of the air column. If we look 
at thickness change as a result of the advection of thickness lines 
only, it appears that there is no heat loss or gain as the air moves 
across the surface of the earth nor are there vertical motions. This 
is not true, of course, since the air gains or loses heat through 
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radiative, dynamical, and conductive processes, and through latent 
heat exchange resulting from phase changes of the water constituent. 

In the model used in this study, in addition to advection, an 
empirical heating term is developed to account for heating or cooling 
by evaporation or condensation and for turbulent transfer of 
sensible heat over ocean surfaces. 

To account for the dynamic processes, vertical velocity is 
assumed to be a second order function of pressure. The 1000-mb 
vertical velocity is composed of a terrain upslope term and a 
friction term. Heat gained or lost by radiative processes was 
ignored as neglibible over the 24-hour period. 

2. Development of the Prognostic Model 

A two layer model is used in the development of the prognostic 
equations. The lower troposphere is divided into two layers as 
illustrated in Fig. 1. 

By use of the equation of state and the hydrostatic equation, 
the thermodynamic equation may be reduced to the following form 




Replacing the vertical derivatives above by finite difference and 
simplifying leads to 
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Where , the thickness. 

In the layer 1000 to 500 mb it has been noted that O) can 
generally be approximated by a second order equation in pressure. 

OJ = 4(x,y,t) -t* u J,o (3) 

By means of equation (2) it may be shown that 



^8 - 4* '/tf ^ 4" UJ, 0 




Taking the vorticity equation at 500 mb in the form 



It +v-v% - 



and substituting from equation (6) gives 



(4) 

(5) 

(6) 



+ v-V^ 5 ) 



Equation (2) is now evaluated for layer 2 giving 



(7) 



tr^Y-Vk + 




The above equation may <be 7 written 



a. 
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Where 



=(£*. [q — ( b i-ri h ^) 

\Cft 6 '^ A f / 

Combing equations (5) and (8) and simplifying leads to 




With the geostrophic assumption for relative vorticity 




it follows that 

= -y + MffO, uo> 

= - V ' n +■ (id 

Subtracting equation (.11) from (10) yields 

% + %(§£),- ^(J $ 0 us) 

However, 

= ( ul>b-UM (W 8 - UJ b ) = -1 (a + b) 

©),„= = * 

thus equation (12) becomes 

(<kfcb) U3) 

dt. is Tff 
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Substituting equations (7) and (9) into equation ( 13 ) and expressing 
advective terms in Jacobian form we obtain the following: 



V a* + ^g3M) 



(14) 



where 



it, 

I 



%Z (Z. I0 ,Z 7 ) = % J (RoA) ;Vo*7V-7(D/4o) 

and 1 

f= to) , g = ( v^n») 

' V Qjl 



The prognostic thickness equation is thus determined fcr the lower 
layer. This equation is of the Helmholtz type, the solution of 
which may readily be obtained by use of a digital computer. 

We' 'Will ho 7 w v . obtain a prognostic equation for the thickness 
of the upper layer. 

Evaluating equation ( 2 ) for the upper level centered at 625 mb 
results in the following: 




csru * + ‘^&6' 



(15) 



where 



[\ 

\c f f b n ‘ 




Substituting from equations ( 4 ), (7) and (9 ) and simplifying, we 



obtain 



^ - %y B -vw - £ 



+IS RAf 
73V fa 



7^ 

cL - 




( 16 ) 
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Substitution of equation (16) into equation ( 15 ), and expressing 
in Jacobian form provides the thickness tendency (prognostic) 
equation for the upper layer of the model. 




( 17 ) 



where 



V-Vh, = %/ 4 T(&,A) 



III. Development of the heating term. 

Initially it will be assumed that the heating rate is expressed 
by an equation of the following format 



Over land the heating term is quite small, and we shall 
consider it as negligible. Over water, excluding radiation, the 
main contributions to the heating term are due to sensible and 
latent heat exchange between the sea and the atmosphere . It is 
of course realized that latent heat gained through evaporation may 
not be utilized by the atmosphere at the exact instant and position 
that it is obtained. However, in the relationship derived, an 
empirical adjustment constant is available should the values 
calculated be determined to be somewhat too small or large. 




( 18 ) 
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According to Laevastu C 3 Zl the heat flux from the ocean to 
the atmosphere due to sensible heat transfer and evaporation or 



condensation may be expressed as follows: 

Qm -M(oAfc+anv,)(T s -TO <w> 

* 

Qfe = {o.2\c+o.i7V s )(. e ises-eo~) ( 20 ) 

^ ft -ToJ) ^ 0 

Qph = ~SV < (Ts *— T <0 ( 21 ) 

Q fc - .<mw(e 5 -CO L (22) 

for (Ts - I a.j'C O 



It is apparent that the contribution from the latter equation 
is small and hence will be considered zero. 

By equating the relationship ( 19 ) given by Laevastu for 
sensible heat flux with the integral of equation (L 8 ) with respect 
to mass per unit area, we may determine a reasonable value for 
the constant "A" . 

Integrating (18) with respect to mass per unit area gives 



U 




where 

Introducing the hydrostatic relationship, p d£- and integrating 

r 

with respect to pressure, and then equating this term to the heat 
flux (19) yields 
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Solving this equation for A, and placing ^U4)- (A^4+0, HV S ) 
leads to a value of 63.5 for A. The resulting term for the sensible 
heat contribution becomes 



Furthermore, since the wind used in the prognostic model will 
be the geostrophic wind at 1000 mb which is about 130 to 140 per- 
cent greater than the surface wind, the constant 0.77 in equation ( 23 ) 
was reduced to 0.5. This constant may be refined to further enhance 
the accuracy of the forecast. 

The final form for this term thus becomes 



Bowen (1926) related the convective transfer of sensible heat 
to evaporation by the following relationship 

Q« = S = K (Ti-Te.) W 

(e s -eO 1000 

A. H. Gordon \_2~\ (1952) using meteorological data from 

British ocean weather ships in the North Atlantic related the mean 




(23) 




( 24 ) 



for (u-u)^o 
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values of vapor pressure differences (e g -e a ) to individual whole 
degree temperature differences of (Ts-Ta) for wind forces 0-3, 4-6 
and greater than 6. 

Roll C4J using this data computed the Bowen ratio from 
corresponding average temperatures and water vapor differences 
(Ta-Ts), (e s -e a ) for wind forces Beaufort 4 and 8. Only a very 
slight increase in the Bowen ratio was observed for force 8 winds. 
The results of those computations are presented in Fig. 2. 

The curve of Fig. 2 was statistically fitted using a linear 
rational transformation and solving for the Bowen ratio in terms of 
temperature difference (Ts-Ta) with the following result: 



P, - Q* — , <^. CT$ -Ta- V 
w ( + . ia.cn -TO 



( 26 ) 



As was previously indicated, the diabatic heating term is to 

• I 

be comprised only of the two elements 0 H and Q E , hence the 

use of relationship (26) enables us to eliminate the vapor pressure 

0 

differences as a parameter in the forecast. Thus Q (Total Diabatic 
Heating Rate) = Q* "+■ Qe 



= I ~t~ . 3 A (T& — ~Te^ ) 



(27) 



The total heating rate determined from equations (2i+) and (27) 
is thus 

Q =3H, S (l 



(28) 
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From equations (2) and (18) , the rate of change of thickness 



due to diabatic heating in the lower and upper layers, which will 
be termed Q 2 and Q-^ respectively, may be expressed as follows 



n = 

^ Vffe 



(ftA4(V s )(T s -Tj 



- 1 !jX !\iU) (t s -t. ) 

Q,= 

Q, = s h Q, 



We have thus established a simplified distribution of the 
effects of the diabatic heating. 

From equations (28) and (31) the rate of change of thickness 
due to diabatic heating in the lower layer is then 

Q a = 7 a a 

= \S.3l (t + .3a(T 5 -T<0)Ui‘>+.S V<jS|o) 

for (r s -TJ?o 



Similarly, it may be shown that 

Q = 0.^1 (v v ,o)(t s -t«J) 

for (J s -To . ) < O 



(29) 

(30) 
(3D 



(32) 

(33) 



20 



In comparison with equation (32) , it is apparent that the 
magnitude of the thickness change due to this term will be quite 
small. However, the constant .27 may be increased if this 
contribution is deemed too small. 

IV. Development of the terrain pressure and friction term 

According to Berkofski and Bertoni C O / the vertical velocity 
produced by terrain may be expressed as follows: 



Making the geostrophic assumption, and also assuming that 
the 1000-mb wind is valid for all levels of terrain, this relationship 
may be expressed as follows: 



It is apparent that this term may be somewhat small in cases 
of high terrain due to the fact that the geostrophic wind at 1000 mb 
is used rather than an integrated wind over the entire height of the 
terrain. However, in nature, the wind will normally tend partially 
to parallel the terrain contours, hence this reduction will somewhat 
compensate for this deficit. 

It was originally intended to use a more refined function for 
terrain induced vertical velocities and to include the effects of 
friction. Time limitations however precluded any further work on 
this term. 




(34) 




(35) 
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It is felt that further expansion of this term or the inclusion 



of the terrain term now calculated by FNWF into the program will 
render some improvement to the prognosis. 

V. Programming the Prognostic Equations for the Control 

Data 1604 Computer 

The prognostic equations were put into finite-difference form 
and programmed for the Control Data 1604 computer used in 
conjunction with the FE880 Univac storage drum. The 1604 computer 
is a stored program, general purpose digital computer with a storage 
capacity of 32,768 words in two indpendent storage units. The 
FE880 storage drum is connected to the 1604 computer by a high 
speed channel which conveys information at a rate of 3,000,000 
bits per second. The drum capacity is 256 blocks of 5,120 
48 bit words per block. 

The program was written in symbolic coded relocatable 
assembly program (SCRAP) for a 63 by 63 hemispheric grid, with 
one hourly time steps. 

An option was written into the program to allow for equating 
the vertical velocity. due to terrain and the heating term to zero 
in order to ascertain their respective effects on the prognosis. 

The 750-mb layer is not considered a standard layer in the 
atmosphere, so no synoptic information is available for the 
initialization of the program. To account for this, the following 
equation developed from the thickness equation was derived. 
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(36) 



D 1S0 = D 100 + h svn5 - [sd.18 (t 8S) - r vo ) + a oil, J 



where h gt775 is the standard thickness of the layers 700 to 750 mb. 

Once the initialization of the program is completed the 750 mb 
height may be computed using the forecast thickness for the upper 
level and FNWF's 500-mb barotropic prognosis. 

The absolute vorticities (' 1 j 5 ( "^ 7 are computed with FNWF's 

subroutine SAR. This routine computes the absolute vorticity over 
a 63 by 63 grid using the following relationships 



\ = a, Sine 


+ 0 U (A,+A^-tA3 +A+- 4 Ao) 


l 1 


Sme 


Where a x = .0746711 




a 2 = .2704348 





- 1 7 

To use this routine, all data must be scaled 2 . The 

absolute vorticity is equated to zero for all points south of 10 N 
latitude . 

After computing the absolute vorticities and reading in the 
various D fields, the Jacobians may be computed utilizing Fleet 
Numerical Weather Facility standard subroutine SAB. The 
subroutine computes the Jacobian in the following manner 

A = jtB.c; 

= (A* £>)(A v c) - (A^B) (A*c) 

The program inserts zeroes in field A at the boundaries. 
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The next term computed is the V ADg . This term is 
computed by Fleet Numerical Weather Facility's subroutine SAD, 
which computes as follows: 

A = 



= 8 : 






This program also zeroes the field at the boundaries. 

To compute the heating term, Q2 we must obtain an hourly 
surface temperature based on the forecast thickness and also 
compute a 1000-mb geostrophic wind. The sea surface temperature 
was assumed to remain constant throughout the period of the 
forecast. 

The equation developed to determine an hourly surface 
temperature as derived from the predicted thickness fields is as 
follows: 

X-c =3‘HS"x I 0 4 [(| 0 .» h».-133k>jU tei^Q^+a.'*7l^) 



After Q ^ is calculated , each point in the grid is checked to 
determine whether it is a land or sea point. At land points Q ? 
is equated to zero; at sea points, as calculated. 

The terrain term U^q is computed using the Jacobian subroutine 
and also checked against the land/sea table. In this case, sea 
points are zeroed and land points left unchanged. 



(37) 
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The terms now calculated are collected and make up the forcing 
function for the Helmholtz equation (Ilf). This equation is solved 
using FNWF's subroutine SAE, providing A for the first hour of 
the prognosis . 

__x 

To solve the Helmholtz equation, VA -BA - C , the program scans 
through field A laterally from the lower left point to the upper right 
point using the extrapolated Liebmann method of computing a new 
point value. The boundary conditions have already been determined 
and stored in the A field, the guess field, which is originally set 
to zero. 

A: : = A, : + X F (AB)c,j — CuS"|^ 

7 L a-* fc-'r+i J 

Rj.4 = A - (A6k.i-Ct.il' 

11 *L a-’-St.i + i J 

When for all points 1*^ € where £ is the criterion for 
convergence or the point where the residual R is considered 
liquidated, the field is assumed relaxed and a proper solution has 
been obtained at each grid point. 

When A has been determined, all the parameters required 
for the computation of Ah^ are available for use in equation (17) . 

The change in height at 1000 mb is then determined by 
subtraction of A hj plus A hg from the 500-mb height change as 
computed by the FNWF barotropic model. 
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6. Results and Conclusions 



Twenty four hour thickness forecasts were produced from 
initial data on 20, 23 and 29 May 1966. The initial program 
contained a number of smoothing passes at many stages of the 
prediction procedure. As a result the initial thickness prognosis 
was very smooth, almost zonal, as shown in Figure 4* Compar- 
ison of this prognosis with the initial and verifying maps, 

Figures 3 and 6, clearly shows that oversmoothing was a major 
source of error. In figures 5, 10 and 13 for the 20th, 23rd 
and 29th of May the results were slightly more encouraging. 

Here smoothing has been reduced somewhat and lower limits 
were placed on absolute vorticity and stability. Still further 
elimination of smoothing throughout the prediction procedure 
gave the thickness and 1000-mb prognosis for the 29 th in 
Figures 15 and 19, which show genuine promise as a useable 
forecast. No 1000-mb analysis existed for the 21st and 24th. 

Generally, the pressure systems moved realistically and 
the final positions were good. However, some excess move- 
ment was observed which might be partially attributed to the 
smoothing still remaining in the model. Some overdevelopement 
of both high and low pressure systems in Figure 15 was also 
observed. However FNWF has a filtering subroutine which in 
addition to filtering out the ultra short waves reduces the 
amplitude of these centers. This subroutine diminishes the 
effect of simply spreading the patterns to adjacent grid points 
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■which is inherent in smoothing subroutines. Smoothing can also 
give rise to fictitious movement of pressure systems. 

By replacing the smoothing with a filtering pass some of 
the difficulties of excess movement and overdevelopement could 
be overcome. 

Another problem is the tendency of some thermal systems, 
known for their persistance, to move with the general circulation, 
an example is the heat low in Southwestern United States. 
Periodic comparison with climatology could be fruitful in reducing 
this tendency. 

Due to a lack of available computer time little consideration 
could be given to the diabatic heating and terrain terms. Further 
refinement and inclusion of these terms could result in some 
improvement of the prognosis. 

The program running time on the 1604 computer for a 24- 
hour forecast was about 30 minutes, including hourly writing 
on magnetic tape of a number of fields which would be unnecessary 
in an operational model. 

The FNWF thickness forecast in Figure 16, which is based 
on simple horizontal advection and numerous constraints, is 
obviously superior to the prognosis resulting from the model 
under investigation... However, the superiority is a result of 
several years of refinement and empirical adjustment. Time 
did not permit such thorough investigation of this research model. 

In conclusion, it has been determined that the assumptions 
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used in the construction of this model are reasonable. With 
proper adjustment and refinement this model should yield a good 
forecast of thickness. This, in turn, should produce a good 
surface prognosis. 
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Layer 1 



Layer 2 



Figure 1. Two Layer division of the atmosphere between 
1000-mb and 500 -mb. h^ and hg are the thickness of the 
layers p-j^ to p^ and p^Q to py respectively. 
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s 

Figure 2. Bowen Ratio vs. Sea-Air Temperature Difference 
(dashed curve) 

Linear rational transformation for positive values. 
Equation (26) (Solid curve). 
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